#include "AddOns/Analysis/Observables/Primitive_Observable_Base.H"

#include "ATOOLS/Org/MyStrStream.H"
#include <iomanip>

using namespace ATOOLS;

namespace ANALYSIS {

  class SFour_Particle_Observable_Base: public Primitive_Observable_Base {  
  protected:

    ATOOLS::Flavour m_flav[4];
    size_t          m_item[4];

  public:

    SFour_Particle_Observable_Base
    (const ATOOLS::Flavour flav[4],const size_t item[4],
     const int type,const double min,const double max,const int bins,
     const std::string &inlist,const std::string &name);
    
    void Evaluate(const ATOOLS::Particle_List &particlelist,
		  double weight=1.,double ncount=1);
    
    virtual bool Evaluate(const Particle *p1,const Particle *p2,
			  const Particle *p3,const Particle *p4,
			  double weight=1.,double ncount=1) const = 0;

  };// end of class SFour_Particle_Observable_Base

  class Four_Particle_Plane_Angle: public SFour_Particle_Observable_Base {  
  public:

    Four_Particle_Plane_Angle
    (const ATOOLS::Flavour flav[4],const size_t item[4],
     const int type,const double min,const double max,const int bins,
     const std::string &inlist);
    
    virtual bool Evaluate(const Particle *p1,const Particle *p2,
			  const Particle *p3,const Particle *p4,
			  double weight=1.,double ncount=1) const;

    Primitive_Observable_Base *Copy() const;
    
  };// end of class Four_Particle_Plane_Angle

}// end of namespace ANALYSIS

using namespace ANALYSIS;

template <class Class>
Primitive_Observable_Base *
GetFourParticleSelector(const Analysis_Key& key)
{
  ATOOLS::Scoped_Settings s{ key.m_settings };
  const auto parameters = s.SetDefault<std::string>({}).GetVector<std::string>();
  if (parameters.size() < 13)
    THROW(missing_input, "Missing parameter values.");
  ATOOLS::Flavour flav[4];
  size_t item[4];
  for (size_t i{ 0 }; i < 4; ++i) {
    const auto kf = s.Interprete<int>(parameters[2 * i]);
    flav[i] = ATOOLS::Flavour((kf_code)std::abs(kf));
    if (kf < 0)
      flav[i] = flav[i].Bar();
    item[i] = s.Interprete<size_t>(parameters[2 * i + 1]);
  }
  return new Class(flav,item,
                   HistogramType(parameters[11]),
                   s.Interprete<double>(parameters[8]),
                   s.Interprete<double>(parameters[9]),
                   s.Interprete<int>(parameters[10]),
                   parameters[12]);
}

#define DEFINE_FOUR_OBSERVABLE_GETTER_METHOD(CLASS,NAME)		\
  Primitive_Observable_Base *					\
  ATOOLS::Getter<Primitive_Observable_Base,Analysis_Key,CLASS>::operator()(const Analysis_Key& key) const \
  { return GetFourParticleSelector<CLASS>(key); }

#define DEFINE_FOUR_OBSERVABLE_PRINT_METHOD(NAME)		\
  void ATOOLS::Getter<Primitive_Observable_Base,Analysis_Key,NAME>::PrintInfo(std::ostream &str,const size_t width) const \
  { str<<"[flav1, item1, ... flav4, item4, min, max, bins, Lin|LinErr|Log|LogErr, list]"; }

#define DEFINE_FOUR_OBSERVABLE_GETTER(CLASS,NAME,TAG)		\
  DECLARE_GETTER(CLASS,TAG,Primitive_Observable_Base,Analysis_Key);	\
  DEFINE_FOUR_OBSERVABLE_GETTER_METHOD(CLASS,NAME)		\
  DEFINE_FOUR_OBSERVABLE_PRINT_METHOD(CLASS)

#include "AddOns/Analysis/Main/Primitive_Analysis.H"

SFour_Particle_Observable_Base::SFour_Particle_Observable_Base    
(const ATOOLS::Flavour flav[4],const size_t item[4],
 const int type,const double min,const double max,const int bins,
 const std::string &inlist,const std::string &name):
  Primitive_Observable_Base(type,min,max,bins)
{
  m_name=name;
  for (size_t i(0);i<4;++i) {
    m_flav[i]=flav[i];
    m_item[i]=item[i];
    m_name+="_"+ToString(m_flav[i])+"-"+ToString(m_item[i]);
  }
  m_name+=".dat";
  m_xmin=min;
  m_xmax=max;
  m_listname=inlist;
}

void SFour_Particle_Observable_Base::Evaluate
(const ATOOLS::Particle_List &inlist,double weight,double ncount)
{
  int no(-1);
  size_t pos[4]={std::string::npos,std::string::npos,
		 std::string::npos,std::string::npos};
  for (size_t k(0);k<4;++k) {
    no=-1;
    for (size_t i(0);i<inlist.size();++i) {
      if (inlist[i]->Flav()==m_flav[k] || 
	  m_flav[k].Kfcode()==kf_none) {
	++no;
	if (no==(int)m_item[k]) {
	  pos[k]=i;
	  break;
	}
      }
    }
  }
  for (size_t k(0);k<4;++k) if (pos[k]==std::string::npos) return;
  Evaluate(inlist[pos[0]],inlist[pos[1]],
	   inlist[pos[2]],inlist[pos[3]],weight,ncount);
}

DEFINE_FOUR_OBSERVABLE_GETTER(Four_Particle_Plane_Angle,
			      Four_Particle_Plane_Angle_Getter,"FourPlaneAngle")
  
Four_Particle_Plane_Angle::Four_Particle_Plane_Angle
(const ATOOLS::Flavour flav[4],const size_t item[4],
 const int type,const double min,const double max,const int bins,
 const std::string &inlist):
  SFour_Particle_Observable_Base(flav,item,type,min,max,bins,
				 inlist,"FourPlaneAngle") {}

bool Four_Particle_Plane_Angle::Evaluate(const Particle *p1,const Particle *p2,
					 const Particle *p3,const Particle *p4,
					 double weight,double ncount) const
{
  Vec4D n1(1.0,cross((Vec3D)p1->Momentum(),(Vec3D)p2->Momentum()));
  Vec4D n2(1.0,cross((Vec3D)p3->Momentum(),(Vec3D)p4->Momentum()));
  p_histo->Insert(n1.Theta(n2),weight,ncount);
  return true;
}

Primitive_Observable_Base *Four_Particle_Plane_Angle::Copy() const
{
  return new Four_Particle_Plane_Angle(m_flav,m_item,m_type,m_xmin,m_xmax,
				       m_nbins,m_listname);
}

